#get markov
get_markov <- function(wd){
  setwd(wd)
  list_file <- list.files()
  list_markov <- list()
  name_markov <- c()
  for (i in c(1:length(list_file))){
    list_markov[[i]] <- read.csv(list_file[i], header = TRUE, sep = ",", quote = "",row.names("cycle"),
                                 dec = ".", fill = TRUE, comment.char = "")
    list_markov[[i]] <- subset(list_markov[[i]],select = c(2:39))
    name_markov[i] <- substr(list_file[i],10,(nchar(list_file[i])-4))
  }
  names(list_markov) <- as.numeric(gsub("[^0-9]", "", name_markov))
  return(list_markov)
}

#get mean markov
get_mean_markov <- function(markov){
  markov <- markov[,c(2:115)]
  col_names<- colnames(markov)
  #extract mean colum
  markov_mean <- markov
  for (i in c(1:length(colnames(markov_mean)))) {
    name_col <- col_names[i]
    if (grepl("LB",col_names[i]) == TRUE) markov_mean[,col_names[i]] <- NULL
    if (grepl("UB",col_names[i]) == TRUE) markov_mean[,col_names[i]] <- NULL
  }
  #change column names
  col_names_mean<- colnames(markov_mean)
  for(i in c(1:length(col_names_mean))){
    col_names_mean[i] <- substr(col_names_mean[i],1,(nchar(col_names_mean[i])-5))
  }
  colnames(markov_mean) <- col_names_mean
  return(markov_mean)
}

#get markov 2
get_markov_2 <- function(wd){
  setwd(wd)
  list_file <- list.files()
  list_markov <- list()
  name_markov <- c()
  for (i in c(1:length(list_file))){
    list_markov[[i]] <- read.csv(list_file[i], header = TRUE, sep = ",", quote = "",row.names("cycle"),
                                 dec = ".", fill = TRUE, comment.char = "")
    list_markov[[i]] <- get_mean_markov(list_markov[[i]])
    name_markov[i] <- i
  }
  names(list_markov) <- as.numeric(gsub("[^0-9]", "", name_markov))
  return(list_markov)
}


#Get YLL
get_YLL <- function(markov,LE,age_init){
  LE_lost <- c((LE-age_init):1)
  death_by_year <- markov$Death
  death_by_cycle_cum <- death_by_year[1:(LE-age_init)]
  death_by_cycle <- rep(x = 0,length(LE_lost))
  for (i in c(2:length(death_by_cycle))){
    death_by_cycle[i] <- death_by_cycle_cum[i]-death_by_cycle_cum[i-1]
  }
  PYLL <- LE_lost*death_by_cycle
  YLL <- sum(PYLL)
  return(YLL)
}

#Get YLL discount
get_YLL_disc <- function(markov,LE,age_init,disc=0.03){
  LE_lost <- c((LE-age_init):1)
  #print(length(LE_lost))
  discount <- exp(-disc*c(1:length(LE_lost)-1))
  #print(length(discount))
  LE_lost_disc <- LE_lost*discount
  death_by_year <- markov$Death
  death_by_cycle_cum <- death_by_year[1:(LE-age_init)]
  death_by_cycle <- rep(x = 0,length(LE_lost))
  for (i in c(2:length(death_by_cycle))){
    death_by_cycle[i] <- death_by_cycle_cum[i]-death_by_cycle_cum[i-1]
  }
  #print(length(LE_lost_disc))
  #print(length(death_by_cycle))
  PYLL <- LE_lost_disc*death_by_cycle
  YLL <- sum(PYLL)
  return(YLL)
}

#get YLD
get_YLD <- function(markov,DW_list){
  z <- 0
  sum_markov <- lapply(markov[,c(1:ncol(markov)-1)],sum)
  for (i in c(1:length(DW_list))){
    z <- z + sum_markov[[i]]*DW_list[[i]]
  }
  return(z)
}

#Get YLD discount
get_YLD_disc <- function(markov, DW_list, disc=0.03){
  discount_YLD <- exp(-disc*c(0:(nrow(markov)-1)))
  DW_by_HS <- c(unlist(DW_list),0)
  DW_matrix <- matrix(rep(DW_by_HS,nrow(markov)),nrow = nrow(markov),ncol = ncol(markov),byrow = T)
  DW_matrix_disc <- DW_matrix*discount_YLD
  DW_matrix_disc_markov <- DW_matrix_disc*markov
  sum_all_disc <- lapply(DW_matrix_disc_markov,sum)
  YLD_disc <- sum(unlist(sum_all_disc))
  return(YLD_disc)
}

get_sum_YLD <- function(x,y){
  
  z <- 0
  for (i in c(1:length(y))){
    z <- z + x[[i]]*y[[i]]
  }
  return(z)
}
#Get DALY
get_DALY <-  function(Markov, LE, age_init, DW_list){
  YLL <-  get_YLL(Markov, LE, age_init)
  sum_M <- lapply(Markov[,c(1:ncol(Markov)-1)],sum)
  YLD<-vector("list",length = length(DW_list))
  names(YLD) <- names(DW_list)
  YLD <- get_sum_YLD(sum_M,DW_list)
  DALY <-  YLL + YLD
  return(DALY)
}

#Get Years of life
get_YL <- function(Markov){
  Markov$Death <- NULL
  YL <- sum(Markov)
  return(YL)
}

#Get Years of life discounted
get_YL_disc <- function(markov, disc=0.03){
  discount_YL <- exp(-disc*c(0:(nrow(markov)-1)))
  alive <- markov
  alive$Death <- NULL
  YL_row <- apply(alive,1,FUN = sum)
  YL_row_disc <- YL_row*discount_YL
  YL_disc <- sum(YL_row_disc)
  return(YL_disc)
}
#Get Years of life saved
get_YLS <- function(Markov_BC, Markov_strat){
  YL_BC <- get_YL(Markov_BC)
  YL_Strat <- get_YL(Markov_strat)
  YLS <- YL_Strat -YL_BC
  return(YLS)
}

#Get Years of life saved, discounted
get_YLS_disc <- function(Markov_BC, Markov_strat,discount=0.03){
  YL_BC <- get_YL_disc(Markov_BC,disc=discount)
  YL_Strat <- get_YL_disc(Markov_strat,disc=discount)
  YLS_disc <- YL_Strat -YL_BC
  return(YLS_disc)
}

get_incidence_CC <- function(dataframe){
  alive <- 0
  alive <- dataframe
  alive <- alive[-1,]
  alive$Death <- NULL
  nb_cirr <- sum(alive$CC.ntrt) + sum(alive$CC.trt)
  alive$CC.ntrt <- NULL
  alive$CC.trt <- NULL
  inc_cirr <- nb_cirr/sum(alive)
  return(inc_cirr)
}

get_incidence_DCC <- function(dataframe){
  alive <- 0
  alive <- dataframe
  alive <- alive[-1,]
  alive$Death <- NULL
  nb_DCC <- sum(alive$DCC.ntrt) + sum(alive$DCC.trt)
  alive$CC.ntrt <- NULL
  alive$CC.trt <- NULL
  inc_DCC <- nb_DCC/sum(alive)
  return(inc_DCC)
}

get_incidence_HCC <- function(dataframe){
  alive <- 0
  alive <- dataframe
  alive <- alive[-1,]
  alive$Death <- NULL
  nb_HCC <- sum(alive$HCC)
  alive$HCC <- NULL
  inc_HCC <- nb_HCC/sum(alive)
  return(inc_HCC)
}